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Many small organisms self-propel in viscous fluids using travelling wave-like deformation of their 
bodies or appendages. Examples include small nematodes moving through soil using whole-body 
undulations or spermatozoa swimming through mucus using flagellar waves. When self-propulsion 
occurs in a non-Newtonian fluid, one fundamental question is whether locomotion will occur faster or 
slower than in a Newtonian environment. Here we consider the general problem of swimming using 
small-amplitude periodic waves in a viscoelastic fluid described by the classical Oldroyd-B constitu¬ 
tive relationship. Using Taylor’s swimming sheet model, we show that if all travelling waves move 
in the same direction, the locomotion speed of the organism is systematically decreased. However, 
if we allow waves to travel in two opposite directions, we show that this can lead to enhancement of 
the swimming speed, which is physically interpreted as due to asymmetric viscoelastic damping of 
waves with different frequencies. A change of the swimming direction is also possible. By analysing 
in detail the cases of swimming using two or three travelling waves, we demonstrate that swimming 
can be enhanced in a viscoelastic fluid for all Deborah numbers below a critical value or, for three 
waves or more, only for a finite, non-zero range of Deborah numbers, in which case a finite amount 
of elasticity in the fluid is required to increase the swimming speed. 


I. INTRODUCTION 

A variety of small prokaryotic and eukaryotic organ¬ 
isms exploit viscous forces from a surrounding fluid in 
order to self propel. The low Reynolds number at which 
they swim means there are no inertial effects, and work 
must be constantly expanded by the cells to produce mo¬ 
tion. This is achieved, for example, by the rotation of 
rigid helical appendages [1] or by the propagation of pla¬ 
nar travelling waves along a flexible flagellum [2]. Our 
fundamental understanding of swimming cells has in¬ 
creased dramatically with the advancement of imaging 
techniques and computer power for more realistic numer¬ 
ical simulations [3]. 

The vast majority of work on swimming at low 
Reynolds number has focused on swimmers moving in 
Newtonian fluids. However, in vivo , many self-propelled 
organisms progress through non-Newtonian fluids. Ex¬ 
amples include the motion of cilia in lung mucus [4], ne¬ 
matodes travelling though soil [5], bacteria in their host’s 
tissue [6], and spermatozoa swimming though cervical 
mammalian mucus [7]. An important question is how a 
transition from a Newtonian to a non-Newtonian fluid 
affects the dynamics and kinematics of micro-swimmers. 
In this paper, we use a simplified modelling approach to 
quantify whether non-Newtonian stresses can help the 
micro-swimmers go faster or if they hinder their motion, 
and how this affects their mechanical efficiency. 

Experimental studies have not yet reached a clear con¬ 
sensus on whether viscoelasticity increases or decreases 
swimming velocities. Instead a range of results has been 
reported for different kinematics and rheological proper- 
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ties. Nematodes swimming in concentrated solutions of 
rod-like polymers undergo an increase in swimming speed 
[8]. In that case, the polymers, aligned by the stress 
caused by the nematode, form local nematic structures 
which give rise to shear-thinning and aid the forward 
propulsion of the nematode. In contrast, solutions of 
long flexible polymers with no shear-thinning but strong 
elasticity lead to a decrease of the nematode’s swimming 
speed [9]. An experiment imitating Taylor’s classic swim¬ 
ming sheet [10] in rotational (planar) geometry shows ex¬ 
actly opposite effects with an increased locomotion in a 
Boger (constant-viscosity, elastic) fluid but a decrease in 
a shear-thinning fluid [11]. Recently, the locomotion of 
flexible-tailed swimmers was also shown to be enhanced 
in a Boger fluid [12]. 

Previous theoretical studies addressing motion in com¬ 
plex fluids have considered a variety of kinematics, in¬ 
cluding undulatory motion [13, 14], helical rotation [15], 
squirming [16, 17], three-sphere models [18], and pad- 
dlers [14], Methods that are ineffectual in a Newto¬ 
nian fluid due to reversibility [19], such as flapping [20] 
or solid body rotation [21], can also be exploited in 
a non-Newtonian setting to induce propulsion [22, 23]. 
In the case of locomotion using helical flagella, small- 
amplitude helices always go slower, but for larger ampli¬ 
tudes, a modest increase is possible [24], In this paper, 
we focus on planar wave motion, a situation for which 
there is a wealth of work starting with Taylor’s swim¬ 
ming sheet [10]. In the presence of a surrounding elas¬ 
tic structure, non-Newtonian stresses were shown com¬ 
putationally to lead to faster and more efficient swim¬ 
ming [25]. Numerical simulations also demonstrated that 
for high-amplitude motion, both shear-thinning [14] and 
polymeric Oldroyd-B fluids [26-28] could lead to faster 
locomotion. In particular, using simulations on finite 
swimming sheets it has been shown that front-back stress 
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asymmetry together with swimmer flexibility leads to in¬ 
creased swimming speeds [28]. 

Analytical work on locomotion by waving focuses on 
small-amplitude motion. In the case of isolated swim¬ 
mers, enhanced swimming was predicted theoretically to 
take place in gels [29], Brinkmann fluids [30] and - with 
the addition of elastolrydrodynamic effects viscoelas¬ 
tic fluids [31], but not in inelastic shear-thinning fluids 
[32]. Two nearby swimmers also synchronise faster in an 
elastic fluid than in a Newtonian medium [33]. However, 
in the case of polymeric fluids, asymptotic results pre¬ 
dicted a systematic decrease of the swimming speed for 
all constitutive models, including all Olclroyd-like fluids 
[13] and general linear viscoelastic fluid models [34] in the 
case of prescribed waveform swimming. A decrease also 
takes place in the case of helical small-amplitude motion 
[15, 35]. Provided the prescribed waving amplitude is 
small compared to its wavelength, it appears thus that 
an isolated swimmer is always slowed down by viscoelas¬ 
tic stresses. 

In this paper, we consider mathematically the most 
general problem for planar locomotion using small- 
amplitude waves periodic both in space and in time. 
Specifically, we prescribe the shape deformation as a sum 
of waves travelling with different wavenumbers and fre¬ 
quencies and in different directions and consider the re¬ 
sulting locomotion in a viscoelastic, Oldroyd-B fluid. We 
show that swimming in a non-Newtonian fluid at small 
amplitudes need not always lead to slower swimming 
compared to the Newtonian case, provided the right com¬ 
bination of waves are considered. For swimming enhance¬ 
ment to be observed, different waves need to travel in 
opposite directions, and the enhancement in that case re¬ 
sults from the asymmetric viscoelastic damping of waves 
with different frequencies. A change of the swimming 
direction is also possible. After presenting the general 
derivations, and introducing a sufficient condition for en¬ 
hanced locomotion, we analyse in detail the cases of two 
or three travelling waves. The enhancement in a vis¬ 
coelastic fluid can be obtained for all Deborah numbers 
below a critical value or, in the case of three waves or 
more, only if a finite amount of elasticity is present in 
the fluid. 


II. GENERAL SMALL-AMPLITUDE WAVE IN 
A VISCOELASTIC FLUID 

A. Setup 

Analogous to Taylor’s classic swimming calculation 
[10, 13], an infinite inextensible sheet of negligible thick¬ 
ness is placed in a fluid and undergoes waving motion. 
The waveform of the sheet is prescribed, and results in 
swimming. In the frame of the swimmer the oscillation 
of the vertical position, y(x,t), of the sheet is described 


by 

+oo +oo 

y(x,t) = b J2 E <*n, m e i{mkx - nut) , (1) 

n——oo m =—oo 

where x denotes the coordinate along the average sheet 
axis and t time. In Eq. (1) the modes n = 0 and to = 0 
are omitted as there is no mean deformation in x or in 
time. The fluid is assumed to be located above the sheet 
along the y > 0 direction In Eq. (1), b is the sheet am¬ 
plitude, k the fundamental wavelength and ui the funda¬ 
mental frequency. We allow both positive and negative 
values of the mode number (to, n) in order to include 
waves travelling in both directions along the sheet. The 
order-one complex coefficients a ra , m represent dimension¬ 
less Fourier amplitude of each (to, n) mode and since y 
is real they satisfy a_ nj _ TO = a* m . To simplify notation 
all sums over n and to from — oo to +00 will be denoted 
with a single summation symbol, YP, n m . 

Upon non-dimensionalising x by k _1 and t by w -1 , 
Eq. (1) becomes 

y(x,t) = eY / (*n,me i{mx - nt \ (2) 

n,m 

with a prefactor e = bk defined as the ratio of the sheet 
amplitude to its wavelength. We assume that this ratio is 
small in this paper, e«l, allowing the swimming speed 
to be computed as an asymptotic expansion in e. 

As the sheet is infinite along the z direction we can 
reduce the three-dimensional swimming problem to two 
dimensions. The velocity field is written as u = u x e x + 
u y e y . This allows a streamfunction, ip(x,y,t), to be de¬ 
fined such that u x = dip/dy and u y = —dip/dx, ensuring 
that the flow remains incompressible. 

In order to find the streamfunction we must first con¬ 
sider the boundary conditions imposed on the flow. On 
the waving sheet the no slip boundary condition enforces 
the velocity of the fluid at the sheet location to be the 
same as the velocity of the sheet, so that 

^^\x, y (x,t) =e^ma„, m e' (mi_nf)e I . (3) 

n,m 

Far away from the sheet we expect that the flow will be 
unaffected by the wavemotion. Hence in the frame of the 
swimmer, the far field velocity will be the speed of the 
swimmer, but in the opposite direction. So if the steady 
swimming of the sheet is denoted —Ue x then we have the 
boundary condition 

V Ip |x ,00 ~ (4) 

where the value of U is to be determined. 

B. Constitutive relationship: Oldroyd-B fluid 

The swimmer is self-propelling in a fluid described by 
the Oldroyd-B constitutive relationship, modelling a di¬ 
lute solution of infinitely extensible polymers in a Newto¬ 
nian solute as a homogeneous continuum [36, 37]. In this 
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classical model, the shear viscosity is constant but the 
polymer elasticity affects the flow, giving rise to normal 
stresses. This is a good model for the Boger fluids used in 
many non-Newtonian micro-swimmer experiments [37]. 
Furthermore, from Ref. [13], we expect that to second 
order in e, our asymptotic results will remain valid for a 
large class of constitutive relationships. 

If p denotes the pressure and r the deviatoric stress, 
Cauchy’s equation of mechanical equilibrium in the ab¬ 
sence of inertia is simply written 

Vp = V • r. (5) 

In an Oldroyd-B fluid, the total deviatoric stress, r, 
a combination of stresses from the Newtonian solvent 
r s , and those from the polymers T p . is written as t = 
Ts + Tp. If tj s denotes the solvent viscosity and assuming 
that T p follows a first-order Maxwell constitutive equa¬ 
tion with relaxation time A, elastic modulus G , and poly¬ 
mer viscosity p p = G/ A, the total stress obeys [37] 

v v 

T + A r= rjj + r]s A 7 , ( 6 ) 

where 7 is the shear rate tensor, defined as 7 = Vu + 
Vu T , and 77 = r] s +rj p is the sum of the solvent and poly¬ 
mer viscosities. In Eq. ( 6 ), the upper-convected deriva¬ 
tive defines the rate of change of the tensor A while it 
translates and deforms with the fluid and is written as 

V f) A 

A= —+u-VA-(Vu t -A + A-Vu). (7) 
at 

Upon non-dimensionalising stresses by pui and shear rates 
by oj , Eq. ( 6 ) becomes 

v v 

r + De r= 7 + /3 De 7 , ( 8 ) 

where /3 = Vs/rj < 1, and De = Aw is the Deborah num¬ 
ber that describes the relative importance of viscoelas¬ 
ticity by comparing the relaxation time to the timescale 
on which the fluid is perturbed, given by 1 /w, where to 
is the fundamental waving frequency. 


1. Solution at order e 

The leading-order constitutive equation is linear and 
given by 

n + + ( 12 ) 

This can be reduced into a streamfunction equation by 
taking its divergence, combining with Eq. (5), and taking 
the curl to eliminate the pressure, leaving 

(l + 0DeJ^ VVr = 0. (13) 

The post-transient solution to Eq. (12) is found using 
Fourier notation and solving the biharmonic equation an¬ 
alytically, leading to 

Vfi =Vv-(l + \m\y)e-^ y e^ mx - nt \ (14) 

z ' m 

n,m 

where the first-order boundary conditions, 

VV’llx.o = 5 ~^ ina n,me l(rnx ~ nt) e x , (15a) 

n,m 

and 

VV’iU.oo = 0, (15b) 

are satisfied. Clearly, the first-order solution is the same 
as the Newtonian case, and as expected there is no swim¬ 
ming at this order. 

2. Solution at order e 2 

At order e 2 , the constitutive equation, Eq. (8), is given 

by 


C. Asymptotic solution 

Since we have e <C 1 we seek to find solutions to the 
stress, streamfunction and velocity in terms of perturba¬ 
tive expansion in e, such that 


( 1+D 4)’*-( 1+ '®4h = 

De(Vuf • n + T\ ■ Vui — ui • Vti) 

— /3De(Vuf • 7 + 7 • Vui — u ■ V' 7 1 ). 

(16) 


= eip! + e 2 tp2 + • ■ • , 

(9) 

r = eri + e 2 T 2 + ..., 

( 10 ) 

U = e 2 U2NN + - 

( 11 ) 


The swimming velocity is expected to be quadratic in e, 
and so we focus on the first and second-order solutions 
(the subscript NN is used as a reminder that the final re¬ 
sult for the swimming speed will quantify non-Newtonian 
swimming). 


Using Fourier notation of the form 

A = ^a ("."*) (17) 

n,m 

for any tensor, vector, or scalar, the first-order constitu¬ 
tive equation, Eq. (12), gives access to the Fourier com¬ 
ponent of the first-order stress as 

~(n,m) _ 1 — m/3Der_(n,m) o\ 

Tl “ 1-mDe 71 • 
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As we are interested in the time-averaged swimming, 
it is sufficient to focus on the time-averaged version of 
Eq. (16). We then use Eq. (18) to express the mean of 
Eq. (16) using the Fourier modes of its right-hand-side, 
and obtain 


value can be found using the time-averaged second-order 
boundary condition, 


d(ip 2 ) 


d y 


x,0 


y^nm|a n;m | 2 , 

n,m 


(26) 


/ \ / • \ De(l — /?) 

M-<7 2 >=£t—D^ x 


(Vuj * • 7 + 7 • VuJ - u* • V7) 


(n,m) 


(19) 


With the first-order streamfunction whose Fourier 
component is 


$ n ’"° = a n , m -( 1 + \m\y)e-^ v e irnx , 
m 


( 20 ) 


we obtain the Fourier modes of the flow velocity, 


u 


(n,m) 


— Qn,.m, & 


H„-M v e im -x 


\m\ 2 y 


— (1 + \m\y)im J ’ 


( 21 ) 


the velocity gradient, 


leading to 



Rewriting Eq. (27) with sums in n and m running from 
1 to 00 only, and using that a- n - m = a* m , leads to a 
simplified expression for the final result as 

U 2 nn — - 'y ' 'y ) nm{ |o nm |<an,—m| ) 

n> 1 m> 1 

( 28) 

where the opposite-sign contributions of waves travelling 
in the +x and —x direction are apparent. 


/1 + fdn 2 De 2 \ 

V 1 + n 2 De 2 ) ’ 


III. A SUFFICIENT CONDITION FOR 
ENHANCED SWIMMING 


= a nrn -e~^ y e imx x 

( —im\m\ 2 y m 2 (l + \m\y)\ , . 

y \m\ 3 y — \m\ 2 im\m\ 2 y J ’ ^ ' 


and the shear stress tensor, 


r (n,m) _ ^-\m\ y Amx ( ~2im\m\ 2 y 2\m\ a y 
^ 1 n,m m \ 2|m| 3 y 2im\m\ 2 y 

(23) 

The divergence and curl are then taken, as before, to 
obtain an explicit equation for the second-order stream- 
function as 


2) y- 1 12 » 2 (P- 1 ) Dc ., 

dy 4 ^ n ' m m 2 1 — inDe 

n,m 

d 2 r 1 

—j 77 e~ 2 ^ y (—4im\m\ 4 y + 4im\m\ 5 y 2 — 2\rri\ 3 im) . 
d y z L J 

(24) 


Integrating with respect to y three times, this gives 


d(V> 2 ) 

dy 


= Ay 2 + By+ C+ 
n 2 (/3 — l)De 


in 2 1 — inDe 


e 2 l m l^(— 2im\m\ 4 y 2 + im|m| 2 ). 

(25) 


Given the form of the boundary conditions at infinity, 
Eq. (4), we obtain A = B = 0 and C is equal to the 
second-order swimming speed, hence C = U' 2 nn- Its 


The result in Eq. (28) gives the leading-order swim¬ 
ming speed of the swimming sheet with the most general 
shape deformation periodic in both x and t. When there 
are no viscoelastic effects De = 0, and the Newtonian 
result is recovered. We denote the swimming speed U 2 N 
in that case. 

As can be seen in Eq. (28), it is the value of the (di¬ 
mensionless) frequency n that affects the non-Newtonian 
change of each mode, not the value of the (dimensionless) 
wavenumber m. In order to gain insight into the condi¬ 
tions for swimming to be enhanced or slowed down by 
the presence of viscoelastic stresses, let us focus on the 
simple case where only the modes \m\ = \n\ are present. 
The sheet deformation is written now as a linear super¬ 
position of travelling waves 

y = e a +n e in ^ + a. n e in ^ , (29) 

ri> 1 


where a+ n and a_ n describes the ?rth mode wave trav¬ 
elling to the right (a: > 0) and left (x < 0) respectively. 
Using Eq. (28) this leads to non-Newtonian swimming 
with speed 


U 2 NN — 1 


n> 1 


/1 + n 2 /3 De 2 \ 

V 1 + n 2 De 2 J ’ 


(30) 


and Newtonian swimming with speed, 

U 2 N = a n , (31) 

n> 1 


where, we have further simplified notation such that 


a n = 2n 2 (\a +n \ 2 -\a- n \ 2 ), (32) 
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describes the superposition of mode n waves in both 
directions. Clearly, for both Newtonian and non- 
Newtonian cases, the addition of backwards waves al¬ 
ways reduces the absolute value of the swimming speed. 
Let us now focus on the relative change in speed when 
comparing swimming between a Newtonian and a non- 
Newtonian fluid. 

Using only inspection we cannot, a priori, define 
a range of Deborah number where we expect to see 
an increase in speed from the Newtonian to the non- 
Newtonian swimming (i.e. U 2 NN/U 2 N > !)• In order to 
look for further insight, we consider the infinite and zero 
Deborah number limits. At zero Deborah number, where 
there are no elastic effects, the ratio of swimming speeds 
is equal to 1. In the limit of large Deborah numbers 
De 1, where elastic effects dominate, it is straightfor¬ 
ward to get from Eq. (30) that U 2 NN/U 2 N = P < 1, and 
thus swimming is always eventually decreased. As the 
value of De increases from zero to infinity, the speed ratio 
could monotonically decrease from 1 to /3, in which case 
no enhancement would be seen, or non-monotonically, 
where enhancement could take place. 

Our numerical simulations indicate that in the cases 
where the speed ratio does go above 1 , then in most cases 
it is always increasing in the neighbourhood of De = 0 be¬ 
fore monotonically decreasing to /? (see numerical results 
in Fig. 1 and discussion below). In order to characterise 
the behavior around De = 0, we can compute deriva¬ 
tives and Taylor-expand the ratio of swimming speeds. 
The first derivative dU 2 NN/dDe evaluated at De = 0 is 
zero because the swimming speed depends quadratically 
on the Deborah number. However, the second derivative 
(the curvature) is non-zero, and is given by 


9 2 U2NN 

<9De 2 


De—0 


^ 2n 2 a„(/3 - 1). 

n> 1 


(33) 


need different combinations of amplitudes and frequen¬ 
cies in order to satisfy the condition in Eq. (34). Hence 
a combination of positive and negative a n values are re¬ 
quired. 

Physically, the increase in swimming speed between 
Newtonian and viscoelastic fluids seen here arises from 
the fact that the damping caused by a non-zero Debo¬ 
rah number affects modes with different frequencies dif¬ 
ferently. Specifically, the damping term of the form 
(1 + n 2 f 3Dc 2 )/(1 + n 2 De 2 ) decreases monotonically with 
n. Modes with higher frequencies are therefore damped 
more than those with lower values of n, which provides 
a mechanism for enhanced swimming. 

For illustration, consider two waves travelling in op¬ 
posite directions with the high-frequency (n) wave trav¬ 
elling along the —x direction ( a n < 0 ) and the low- 
frequency (to) one along the +x direction (a m > 0 ). 
Then their respective amplitudes be such that the result¬ 
ing Newtonian swimming speed is positive, U 2 N > 0. In 
the viscoelastic fluid, the a n wave will be damped more 
than the a m wave, as n > to. On one hand, decreasing 
the magnitude of the a n wave will increase the swimming 
speed while on the other hand, decreasing the a m mode 
will hinder the swimming velocity - it is thus a matter of 
relative decrease. If the wave amplitudes are such that 
the gain found by suppressing the a n wave more than 
compensates for the damping of the a m wave, then the 
non-Newtonian swimming speed will be above the New¬ 
tonian one, U 2 NN > U 2 N■ If the wave amplitudes are 
such that U 2 N < 0, then a similar reasoning might be 
used to lead to U 2 NN > 0 and in that case, viscoelastic¬ 
ity might lead to a reversal of the direction of locomotion. 


IV. SUPERPOSITION OF TWO TRAVELLING 
WAVES: CONTINUOUS ENHANCEMENT 


When it is divided by the Newtonian swimming speed, 
Eq. (31), the above gives access to the curvature of 
U 2 NN/U 2 N at De = 0 (this is equivalent to taking the 
first derivative of the speed ratio with respect to Dc 2 ). 
If that curvature is positive, then faster swimming oc¬ 
curs in the neighbourhood of De = 0. As we always have 
(3 < 1, the curvature is positive if there is a sign difference 
between the sums in Eqs. (31)-(33) and therefore a suffi¬ 
cient condition for enhanced swimming is the kinematic 
condition 


Y, a n 

X 

^ n 2 a n 

n> 1 


n> 1 


In order to achieve the condition in Eq. (34), waves trav¬ 
elling in opposite directions are required. Indeed, for 
example if all a n amplitudes are positive, then it is easy 
to see from Eq. (30) that each a n mode decreases in am¬ 
plitude, resulting in an overall decrease in magnitude of 
the speed. If there are waves travelling in both direction, 
i.e. at least one a_ n 7 ^ 0 and one a+ n 7 ^ 0 , then they 


We now consider in detail simple cases. We start by 
swimming using two travelling waves, and show that in 
this case the sufficient condition described above is in 
fact necessary: when enhancement takes place, it will 
lead to faster swimming for all Deborah numbers be¬ 
low a critical value. In order to analytically describe 
situations where faster swimming can occur, two simple 
waveforms each containing two waves travelling in op¬ 
posite directions will be considered. Clearly these two 
travelling waves must have different frequencies, other¬ 
wise they are both damped in the same proportion by 
viscoelasticity and the swimming speed decreases. 

A. Superposition of two travelling waves with 
identical wave speeds 

An example of two waves with different frequencies 
modes, amplitudes, and wave direction but identical 
magnitude of wave speed is given by 

■y(x, t) = e [a sin(x — t) + sin n(x + t )], (35) 
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FIG. 1: Ratio between the non-Newtonian swimming speed U 2 NN, and the Newtonian value U^n, as a function of the Deborah 
number, De, for various values of the relative wave amplitude a, and frequency ratio n, in the waveform from Eq. (35). Here 
we have chosen 0 = 0.1. Left: fixed value of n = 2 and a range of a values (between 0.5 and 5). Right: fixed value of a = 9.5 
and n ranging between n = 1 and 9. 


where a is the dimensionless ratio of amplitudes between 
the two waves. Using Eqs. (30) and (31) for the si¬ 
nusoidal waveform in Eq. (35) we get the second-order 
Newtonian swimming speed as 


U 2 N 



(36) 


while the second-order non-Newtonian swimming speed 
is given by 


rr a 2 f 1 + 0De 2 \ n 2 f 1 + n 2 /3De 2 \ 
2NN ~ Y ^ 1+ De 2 ) ~ Y 1 + n 2 De 2 ) ' 


(37) 


To find where faster swimming occurs, we compute 
as above the second derivative of the swimming ratio, 
U 2 NN/U 2 N with respect to De at De = 0, giving 


d 2 fU2 nn\ 
<9De 2 \ U 2 N ) 


De=0 


2 ( 0 - 1 ) 



(38) 


This is positive (i.e. upwards curving from U 2 NN/U 2 N = 
1) when n < a < n 2 . Hence faster swimming requires the 
relative amplitude between the two waves to lie in a pre¬ 
cise interval. If a is too small the behavior is dominated 
by the —x wave while if it is too large the dynamics is 
dominated by the + 2 : wave. At higher modes, the range 
of amplitudes available to the swimming sheet that can 
produce faster swimming in a non-Newtonian environ¬ 
ment compared to a Newtonian one is increased. 

We illustrate in Fig. 1 these results numerically. We 
plot the ratio of swimming velocities, U 2 NN/U 2 N, as a 
function of the Deborah number, De, for a range of val¬ 
ues of both n and a. We choose a fixed value of 0 = 0.1. 
The computational results confirm that when enhanced 
swimming is obtained, the speed ratio first increases in 
the neighbourhood of De = 0 before monotonically de¬ 
creasing to 0. This validates the curvature analysis as 
a proxy for predicting enhanced swimming, and indeed 
faster swimming in a non-Newtonian fluid is seen in the 


range n < a < n 2 . An illustration of travelling wave that 
swims faster in a non-Newtonian fluid is shown in Fig. 2, 
with n = 2 and a = 5/2. This waveform corresponds to 
the speed ratio shown as the uppermost solid grey line in 
Fig. 1 with a maximum of U 2 NN = 1-3 at De = 0.5. 

Further analytical insight can be provided by noting 
from Fig. 1 that the peak swimming speed ratio occurs 
when De is order one. Dividing the result in Eq. (37) by 
that in Eq. (36) and taking a first derivative with respect 
to De we can compute the value of the Deborah number 
at which the velocity ratio is extremized. It occurs for 
two values of De given by 


/ - Q, 

Dei* = \ ~~ : ri - T\ an d De 2 * = 0. (39) 

U n z {a — 1) 


For a above n 2 the only solution is the maximum value 
of 1 occurring at De 2 » = 0. When a crosses below n 2 a 
maximum is created near Dei* = 0, and increases as a 
decreases. When a = n a transition occurs where Dei* 
changes from a maximum point (n < a) to a minimum 
(a < n)\ its value at that point is Dei* = 1/n. It remains 
a minimum until a crosses the value 1, below which the 
only solution is the maximum of 1 at De 2 * = 0. 

A final point of interest in Fig. 1 is the fact, as discussed 
above, that the ratio between the swimming speeds can 
become negative. In these cases, the swimmer would then 
swim in different directions in the Newtonian and non- 
Newtonian fluids, as was already noted in Ref. [15]. This 
occurs when a < n, and the speed ratio goes through a 
minimum before increasing back towards 0 at large Deb¬ 
orah numbers. The reversal of swimming occurs when 
there is a difference in sign between Eq. (36) and Eq. (37), 
which corresponds to the amplitude range 


n 2 (l + 0n 2 De 2 )(l + De 2 ) 
(l + n 2 De 2 )(l + 0De 2 ) 


(40) 


This result is reminiscent of a recent study on reciprocal 
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FIG. 2: Illustration of a waveform producing faster swimming in a non-Newtonian fluid. The waveform is described by Eq. (35) 
with a = 5/2 and n = 2, and corresponds to a swimming speed ratio as shown by the uppermost solid grey line in Fig. la. 
The black lines in each of the four figures show the waveform at dimensionless times 0, 7r/2, 7r and 37t/2, respectively, and the 
grey lines show the evolution of the wave an eighth of a period later, to show how the wave travels and changes shape during 
its period. 


(time-reversible) motion in a worm-like micellular solu¬ 
tion, which showed that the direction and the speed of 
the swimmer could be changed when distinct Deborah 
numbers are reached [38]. Finally, we can also find a 
range of a values for which the swimmer will not only 
change direction but will also swim with a larger magni¬ 
tude, which occurs when 

< a < n. (41) 

Here the swimming speed ratio becomes negative and less 
than — 1 . 



This condition requires n < a < n 2 , and defines the 
range of Deborah number where forward swimming en¬ 
hancement is achieved, namely [0, De a ]. If we enforce the 
curvature to be negative then we cannot find a set of vi¬ 
able parameters for which U 2 NN > U 2 N > 0 , thus show¬ 
ing that the sufficient condition is also necessary when 
two modes are considered: in the case of two waves, if 
forward swimming enhancement is ever to be obtained, it 
will take place for any Deborah number below a critical 
value De 0 . 


B. Necessary vs. sufficient condition for enhanced 
swimming 

The sufficient condition for enhancement derived in 
§111 detailed the conditions required for an upwards curv¬ 
ing of the swimming speed ratio from zero Deborah num¬ 
ber. In order to study if this sufficient condition is also 
necessary, we search analytically for the conditions lead¬ 
ing to U 2 nn > U 2 n, leading to 



C. Swimming efficiency 

We now turn to energetic considerations. The rate of 
viscous dissipation in the fluid as the sheet is swimming 
is equal to the volume integral of t : 7 in the fluid. 
At leading order we therefore have to integrate T\ : 7 -,. 
With the general waveform in Eq. (1), the dimensional 
second-order dissipation rate in the non-Newtonian fluid 
per unit length in the e z direction is easily found and we 
obtain W = e 2 W 2 NN + ■ ■ ■ with 




FIG. 3: Ratio of the swimming efficiency in a non-Newtonian fluid compared to its Newtonian counterpart as a function of De: 
(a) n = 2 for a range of values of a; (b) a = 9.5 and 1 < n < 9. The waveform is the one described in Eq. (35). 


W 2 nn — EE 8TTrjoj 2 mn 2 

n> 1 ra> 1 


/1 + n 2 /3 De 2 \ 
V 1 + n 2 De 2 / 


(| &n,m 


(43) 


The result in Eq. (43) should then be compared with its 
Newtonian counterpart. 

Let us consider for illustration the waveform in 
Eq. (35). In that case, the ratio of the work done against 
the non-Newtonian fluid compared to the Newtonian one 
is given by 


W 2 nn V 

_ — _ ^ 

W 2 n r/N{n 3 +a 2 ) 

2 ( 1 + Dc 2 77s/77 ^ 3 ( l + n 2 De 2 T] s /r] \ 

a V 1+De 2 J +? V l + n 2 De 2 )\ ’ 

(44) 


We plot the ratio of efficiencies against De for a range 
of relative wave amplitude a, and wavenumber ratio 
n, in Fig. 3, where 773/77 = /3 = 0.1. Clearly, an in¬ 
crease in swimming speed is correlated with an increase 
in efficiency, but increased efficiencies can in fact be ob¬ 
tained without enhanced swimming. Indeed, increased 
efficiency is obtained as soon as 


/ U2NN \ 2 VNW2NN 

V U 2 n ) VW 2 N 


(47) 


Given that the right-hand side of Eq. (47) is less than 
one, the condition for enhanced efficiency does not re¬ 
quire enhanced swimming. Specifically, using the illus¬ 
trative sinusoidal waveform of Eq. (35), we obtain an 
improved efficiency when 


(n 3 + a 2 ) 
( a 2 — n 2 ) 2 


> 


a 


1 + /3De" 
1+De 2 


3^1 + fin 2 De^ 
1 + n 2 De 2 


a 


1 + /3De 2 
1 + De 2 


— n 


1 + f3n 2 De 2 
1 + n 2 De 2 


(48) 


where t]n is the Newtonian viscosity. In order to contrast 
the locomotion in the polymeric fluid with that in the 
solvent alone, we then take tjn = ?y s . Furthermore, as is 
done traditionally, the swimming efficiency is defined as 


rjU 2 

~W' 


(45) 


In order to compare the efficiency of swimming in the 
different fluids, we compute the ratio 


&2NN _ V^2NN W 2 N 
82N W2NN ijnU 2 n 


(46) 


The ratio of the work and viscosity in the two different 
fluids, t]W 2 n/vnW 2 nn, is always greater than 1 for non¬ 
zero Deborah number, meaning that when the swimming 
speed ratio U 2 nn/U 2 n is greater than 1 , the swimming 
efficiency is automatically always increased. 



FIG. 4: Three example waveforms are shown for which the 
swimming speed is not enhanced but the efficiency is. The 
relative amplitude in Eq. (35) lies outside the range n < a < 
n 2 (n = 2 ). 
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FIG. 5: Illustration of a waveform from Eq. (49) with a = 2 and n = 2 that produces faster swimming in a non-Newtonian 
fluid. The black lines in each of the four figures show the waveform at dimensionless times 0, n/2, n and 37r/2, respectively, 
and the grey lines show the evolution of the wave an eighth of a period later. 


for which n < a < n 2 is not a necessary condition. This 
result is illustrated in Fig. 4 in the case n = 2. When a > 
n 2 , the waveform travels in the same direction in both 
fluids and the swimmer is always faster in a Newtonian 
fluid although it is more efficient in the non-Newtonian 
one for a range of Deborah numbers. 


D. Two waves with identical wavelengths 


Instead of two waves with identical wave speeds, en¬ 
hanced swimming can also be obtained in a combination 
of waves with identical wavelengths. Since the waves need 
to have different frequencies, then they necessarily have 
different wave speeds. As an example we consider here 
the waveform 


y = e[a sin(a: — t) + sin(a; + nt)\. (49) 


This gives 


U2NN = vr 


2 V 1 + De 2 ) 2 \1 + n 2 De 


a 2 f 1 + f3De z \ n (1 + n 2 De 2 


, (50) 


and 


C/ 2 at — -(a ? — n). 


(51) 


Similarly as above, the second derivate of U 2 NN/U 2 N at 
De = 0 is given by 


d 2 fU2 nn\ 

DDc 2 \ U‘2N ) 


De—0 


03-1), 


(52) 


and faster swimming occurs when 


n 1 / 2 < a < n 3 / 2 , 


(53) 


which is confirmed by numerical computations (not 
shown). A waveform leading to enhanced swimming in 
this case is illustrated in Fig. 5, in the case a = 2 and 
n = 2. This corresponds to a maximum speed enhance¬ 
ment of U 2 NN/U 2 N ~ 1-1 at Deborah number De « 0.4. 
To obtain the optimal Deborah number, we extremise the 
ratio of swimming speeds to find the peaks occurring at 


Dei* = \ ~Yj —and De 2 * = 0, (54) 

V )r(ayn — 1 ) 


with a behavior qualitatively similar to that of the last 
section. 


V. SUPERPOSITION OF THREE TRAVELLING 
WAVES: CONTINUOUS VS. DISCRETE 
ENHANCEMENT 

In the previous sections, where the superposition of 
two waves was considered, we saw that when enhance¬ 
ment is present, it is continuous from De = 0 to an order 
one Deborah number De 0 , where the value De 0 7 ^ 0 is 
the only non-zero solution of U 2 NNIU 2 N = 1- We now 
demonstrate that if the swimmer is able to use a third 
travelling wave, it is possible for swimming enhancement 
to occur only when a finite amount of viscoelasticity is 
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FIG. 6: Two different types of enhancement are shown for two different three-mode waves with m = 1, ri2 = 4 and 713 = 8: (a) 
continuous enhancement from zero Deborah number in the range [0,De a ] with De a ~ 1.9 ( a ni = 1, a„ 2 = —1, a n3 = 0.3); (b) 
enhancement in a discrete, finite, range of Deborah numbers, [De&,De c ], where Deb ~ 0.3 and De c ~ 1.8 (a ni = 1, a n2 = —2 , 
a n3 — 1 . 2 ). 


present, i.e. for values of the Deborah number in the 
range [Def,,De c ], where Deb and De c are both non-zero. 

We consider, for illustration purposes, a waveform with 
three modes 1 < ri 2 < 713. The corresponding square 
amplitudes ai, a n2 and a n3 are non-dimensionalised by 
ai so that we take ai = 1. Using the same notation as 
above, the Newtonian swimming speed is then given by 

U 2 N = 1 + Ura 2 + a n 3 - (55) 

The difference between the non-Newtonian and Newto¬ 
nian swimming speeds is then found to be 


Uznn — U2N =De“(/3 — 1) 


1 


1 + De 2 


( nla n2 \ 
V1 + n^De 2 ) 


f nj an 3 \ 

\ 1 + n|De 2 ) 

(56) 


Focusing on cases where U 2 N > 0, enhanced forward 
swimming is found when Eq. (56) is positive. As shown 
in Fig. 6 numerically for 712 = 4 and 713 = 8, there are 
two types of enhancements possible: either on a range 
[0,De a ] where the velocity ratio curves upward at the 
origin (continuous enhancement, Fig. 6a, as in § IV) or 
on a range [De&, De c ] for which the curvature at De = 0 is 
initially negative before curving upward as the viscoelas¬ 
ticity increases (discrete enhancement, Fig. 6b). 

In order to distinguish between them analytically, we 
observe that when the curvature is negative, we can ei¬ 
ther have no enhancement or enhancement at a finite 
Deborah number. Hence, we need to search for cases 
where Eq. (56) is positive, given that the curvature at 
the origin is negative. The curvature of the general wave 
was obtained in Eq. (33), hence for negative curvature in 
our three-mode waveform we require 

k = 2(/3 - 1)(1 + n|a„ 2 + n 2 a„ 3 ) < 0. (57) 


The result in Eq. (56) can then be written in terms of 
U 2 n and k as 


n!niDe 4 (/3 - 1) ( ” 2 — 1 + a„ 2 + > 

\ n 2 n 3 J 

[De 2 «(i + De 2 ) + n\n\ De 4 (/3 - 1 )U 2N ] ■ 

(58) 

As f3 — 1 < 0, and assuming that U 2 N > 0 and k < 0, the 
minimum requirement for non-continuous enhancement 
is 

| 

— - 2^2 -k a fi 2 + a n 3 < 0. (59) 

n 2 n 3 

The three conditions given by Eqs. (55) (U 2 N > 0), 
(57), and (59) can be satisfied simultaneously only when 
a n2 < 0 and a n3 > 0, i.e. the first and third modes must 
have a different sign to the second mode. We then search 
numerically over the domain {a n2 < 0, a„ 3 > 0} and De 
to obtain the regions where Eq. (56) is positive provided 
U 2 N > 0 and k < 0, in the example case 712 = 4 and 
n-3 = 8. The values of a„ 2 and a n3 fitting these condi¬ 
tions are shown in Fig. 7 (grey scale domain) while the re¬ 
gion showing continuous enhancement is shown in black. 
The grey scale colouring scheme used in Fig 7 displays 
the value of the lower bound in the interval, Deb, from 
low (dark grey) to high (light grey) values. For three 
waves, in contrast to the case of two waves, situations 
exist therefore where a finite amount of viscoelasticity 
is required to get enhanced propulsion, De > Deb > 0 
(analysis for four and five mode waves show similar re¬ 
sults and are not shown here). 


VI. DISCUSSION AND CONCLUSION 

Motivated by the non-Newtonian environment in 
which many swimmers propel themselves in vivo , in this 
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FIG. 7: Regions in the parameters space {a„ 2 ,a n 3 } where enhanced swimming occurs. The black section represents where 
the curvature is positive (i.e. upwards curving) so that from infinitesimally small Deborah number we get an increase in the 
swimming speed. In contrast, the grey scale section shows regions where enhanced swimming is obtained despite negative 
curvature at the origin. The grey scale color scheme quantifies the value of non-zero Deborah numbers at which the increase in 
swimming speed first occurs, Deb, from low Deborah number in dark, to high Deborah number (ss 0.8) in light grey. Results 
are shown for 72,2 = 4 and 223 = 8 . 


paper, we have calculated the speed of Taylor’s swimming 
sheet in a Newtonian and an Oldroyd-B (non-Newtonian) 
fluid, in the small-amplitude limit. In contrast to pre¬ 
vious analytical studies, we found that small-amplitude 
travelling waves can produce faster swimming in a non- 
Newtonian fluid compared to a Newtonian fluid when 
there are waves travelling in opposite directions in dif¬ 
ferent frequency modes and with different amplitudes. 
Physically, in a non-Newtonian fluid the waves in higher 
frequency modes are damped more than those in lower 
frequency modes, increasing the overall speed of the wave 
under conditions placed upon the difference in frequency 
and amplitude of the summed waveforms. The efficiency 
of the wave can also be increased, and the direction of 
swimming can sometimes be reversed. By studying in de¬ 
tail the superposition of two or three travelling waves, we 
also showed that the range of Deborah number in which 
the enhancement of the swimming speed takes place can 
either include the origin, in which case any small amount 
of viscoelasticity will lead to faster swimming, or it may 
be a finite interval which does not include the origin, 
meaning that faster locomotion requires a finite amount 
of viscoelasticity. 

The results in this paper are reminiscent of recent 
experimental and theoretical work on the role of iner¬ 
tia in locomotion, where two important questions have 
been addressed: (1) for a non-swimmer at zero Reynolds 
numbers, how much inertia is needed to make it swim? 
and (2) how does the locomotion speed of a Stokesian 
swimmer vary with inertia? The answer to question 
(1) depends crucially on the geometry and actuation 
of the swimmer and both discrete [39, 40] and contin¬ 
uous [41] transition to swimming were obtained. In re¬ 
sponse to question (2), model organisms called squirmers 
were shown to vary their speed monotonically from zero 
Reynolds number [42, 43]. Similarly, in our results, we 


showed that a careful design of the swimming kinematics 
could lead to either a decrease or an increase, which could 
be continuous or discrete, of the locomotion speed. We 
expect that these results will remain valid for more realis¬ 
tic models of swimming organisms, in particular those in¬ 
cluding features such as large-amplitude, finite-size, and 
three-dimensional effects. 

In three-dimensions, the same frequency-dependent 
damping term as the one in Eq. (28) is present for in¬ 
finite cylindrical swimmers [35], hence similar results are 
expected to hold. With regards to finite sized swimmers, 
backward propagating waves are expected to occur due 
to the finite nature of real flagella. Previous computa¬ 
tional studies have shown that the addition of viscoelas¬ 
ticity decreases the backwards motion of a finite swim¬ 
mer [26] due to the presence of a viscoelastic network 
behind the swimmer. The opposite has been observed 
experimentally for nematodes where hyperbolic stresses 
created along the swimmer hinder propulsion [9]. It is 
yet unclear how our theoretical results would extend to a 
finite swimmer, though we may expect that the mecha¬ 
nism provided in this paper would provide an additional 
contribution to the swimming speed. In the case of multi¬ 
ple swimmers, it would be interesting to investigate how 
waves with both high and low frequency modes affect 
one another and potentially synchronise. From Ref [33] 
we expect the synchronisation rate to increase with the 
frequency of the waves however the generalisation to mul¬ 
tiple modes in viscoelastic fluids has yet to be done. A 
recent study address a related issue in the Newtonian 
case [44]. 

The waveforms produced here offer insight into how 
swimming speeds can be increased in fluids with vis¬ 
coelastic properties often found in nature [45]. Can these 
shape kinematics occur in biology? For flagellar swim¬ 
ming this requires understanding of how the stochastic 
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actuation of molecular motors create waveforms. Dynein, 
the motor protein causing flagella bending, has been pro¬ 
posed to have two distinct modes to create oscillatory 
bending - these can be described as active and passive, 
or forward and reverse active modes [46], leading to trav¬ 
elling waves that can propagate up or down the flagellum. 
Due to the finite nature of flagella the wave is reflected 
back off the tail end or basal body, thus creating passive 
backwards waves [47]. Experiments on Drosophila sper¬ 
matozoa show that the cells use actively created forwards 
and backwards flagellar waves to avoid obstacles [48]. 
Furthermore by solving elastohydrodynamic force bal¬ 
ance equations on infinite flagella analytic studies have 
shown two different modes of waves travel along the flag¬ 
ella with the same frequency, but different amplitudes 
and directions [49]. Hence a flagellum naturally creates 
forward and backward travelling waves with different am¬ 
plitudes. The enhancement described in this paper re¬ 
quires however waves with different frequencies travel¬ 
ling in the backwards direction for enhancement to oc¬ 
cur. The addition of higher frequency modes (n = 2 and 
n = 3), found in small amounts in beating spermatozoa 
[50], would not however lead to an increase in swimming 
speed as a n > 0 for all values of n found experimentally. 
Changes in flagella beating frequency can occur by alter¬ 
ing the environment in which the swimmer propagates, 
for example hyperactivation when mammalian sperma¬ 
tozoa reach the ovum, leading to a reduction in the beat 
frequency and increase in the beat amplitude [51]; a vari¬ 
ation in ATP or salt concentrations also change frequen¬ 


cies [46]. Recent work on the unsteady modes of flagellar 
motion show that most modes have a frequency smaller 
than the fundamental frequency, and hence would corre¬ 
spond to a reduced swimming speed [52]. The addition of 
noise to the molecular motor oscillations, either through 
variations in concentrations in the bulk or variations be¬ 
tween motors, could lead to increased swimming provided 
the coherent noise is large enough for the flagellum to ac¬ 
cess a higher frequency mode, however this is larger than 
the noise measured [53]. Backwards travelling wave re¬ 
sults have been described for muscle-actuated planar mo¬ 
tion occurring for example in the nematode Caenorhab- 
ditis elegans [54] as well as other flagellar systems such 
as the green alga Chlamydomonas reinhardtii [55]. 

While our study offers only an idealised view, and al¬ 
though as of yet there are no experimental studies from 
biology for which the model here would predict faster 
swimming, our work address the most general periodic 
waving deformation and points to the use of multiple 
waves travelling in different directions as a mechanism 
allowing control of swimming magnitude and direction 
in complex environments. 
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